library(dynlm)
## Warning: package 'dynlm' was built under R version 4.0.5
## Loading required package: zoo
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
library(zoo)
library(forecast)
## Warning: package 'forecast' was built under R version 4.0.4
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
library(strucchange)
## Warning: package 'strucchange' was built under R version 4.0.5
## Loading required package: sandwich
library(tidyverse)
## -- Attaching packages --------------------------------------- tidyverse 1.3.0 --
## v ggplot2 3.3.2 v purrr 0.3.4
## v tibble 3.0.4 v dplyr 1.0.2
## v tidyr 1.1.2 v stringr 1.4.0
## v readr 1.4.0 v forcats 0.5.0
## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x stringr::boundary() masks strucchange::boundary()
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
library(lubridate)
##
## Attaching package: 'lubridate'
## The following objects are masked from 'package:base':
##
## date, intersect, setdiff, union
library(MLmetrics)
## Warning: package 'MLmetrics' was built under R version 4.0.5
##
## Attaching package: 'MLmetrics'
## The following object is masked from 'package:base':
##
## Recall
library(xts)
##
## Attaching package: 'xts'
## The following objects are masked from 'package:dplyr':
##
## first, last
library(DataCombine)
## Warning: package 'DataCombine' was built under R version 4.0.5
library(TSstudio)
## Warning: package 'TSstudio' was built under R version 4.0.5
cat("\014")
set.seed(123)
#clear global environment
rm(list = ls())
#below function reads the data from csv file
#Input : String of file
#Output: dataframe
read_data <- function(data){
data <- read.csv(data)
#make first column as date
data[,1] <-as.Date(data[,1],format ="%m-%d-%y")
data$M2 <- as.numeric(data$M2)
return(data)
}
df<-read_data("Data.csv")
#convert cpi to inflation and make new dataframe storing in df
df <- change(data = df,Var = "CPI",NewVar = "Inflation",type = "percent")
##
## Remember to put data in time order before running.
##
## Lagging CPI by 1 time units.
#convert all columns to monthly time series
inflation <- ts(df$Inflation,frequency=12, start=c(1991,1))
ipilsm <- ts(df$IPILSM,frequency=12, start=c(1991,1))
tbr <- ts(df$TBR,frequency=12, start=c(1991,1))
m2 <- ts(df$M2,frequency=12, start=c(1991,1))
exports <- ts(df$Exports,frequency=12, start=c(1991,1))
imports <- ts(df$Imports,frequency=12, start=c(1991,1))
#split data for train and test and set the last 12 months as a testing partition
split_inflation <- ts_split(ts.obj = inflation, sample.out = 12)
split_ipilsm <- ts_split(ts.obj = ipilsm, sample.out = 12)
split_tbr <- ts_split(ts.obj = tbr, sample.out = 12)
split_m2 <- ts_split(ts.obj = m2, sample.out = 12)
split_exports <- ts_split(ts.obj = exports, sample.out = 12)
split_imports <- ts_split(ts.obj = imports, sample.out = 12)
#fit auto.arima on each series training data
inflation_fit <- auto.arima(split_inflation[["train"]])
ipilsm_fit <- auto.arima(split_ipilsm[["train"]])
tbr_fit <- auto.arima(split_tbr[["train"]])
m2_fit <- auto.arima(split_m2[["train"]])
exports_fit <- auto.arima(split_exports[["train"]])
imports_fit <- auto.arima(split_imports[["train"]])
#generate predictions for next 12 months
inflation_predictions <- forecast(inflation_fit,h=12,fan = TRUE)
ipilsm_predictions <- forecast(ipilsm_fit,h= 12,fan = TRUE)
tbr_predictions <- forecast(tbr_fit,h=12,fan = TRUE)
m2_predictions <- forecast(m2_fit ,h=12,fan = TRUE)
exports_predictions <- forecast(exports_fit,h=12,fan = TRUE)
imports_predictions <- forecast(imports_fit,h=12,fan = TRUE)
#plot all the forecast with actual data with their fanplots
plot_forecast(inflation_predictions)
plot_forecast(ipilsm_predictions)
plot_forecast(tbr_predictions)
plot_forecast(m2_predictions)
plot_forecast(exports_predictions)
plot_forecast(imports_predictions)
#check residuals of all fitted models and conclude by Ljung-Box test
checkresiduals(inflation_fit)
##
## Ljung-Box test
##
## data: Residuals from ARIMA(2,0,0)(1,0,2)[12] with non-zero mean
## Q* = 27.369, df = 18, p-value = 0.07232
##
## Model df: 6. Total lags used: 24
checkresiduals(ipilsm_fit)
##
## Ljung-Box test
##
## data: Residuals from ARIMA(1,0,2)(0,1,2)[12] with drift
## Q* = 36.181, df = 18, p-value = 0.006689
##
## Model df: 6. Total lags used: 24
checkresiduals(tbr_fit)
##
## Ljung-Box test
##
## data: Residuals from ARIMA(0,1,0)
## Q* = 23.149, df = 24, p-value = 0.511
##
## Model df: 0. Total lags used: 24
checkresiduals(m2_fit)
##
## Ljung-Box test
##
## data: Residuals from ARIMA(2,1,2)(0,1,1)[12]
## Q* = 41.622, df = 19, p-value = 0.001993
##
## Model df: 5. Total lags used: 24
checkresiduals(exports_fit)
##
## Ljung-Box test
##
## data: Residuals from ARIMA(2,1,0)(0,0,2)[12]
## Q* = 35.443, df = 20, p-value = 0.01787
##
## Model df: 4. Total lags used: 24
checkresiduals(imports_fit)
##
## Ljung-Box test
##
## data: Residuals from ARIMA(0,1,1)(0,0,1)[12]
## Q* = 42.508, df = 22, p-value = 0.005421
##
## Model df: 2. Total lags used: 24
#plot training with test data
inflation_fit %>%
forecast(h=12) %>%
autoplot() + autolayer(split_inflation[["test"]])
## Warning: Removed 1 row(s) containing missing values (geom_path).
ipilsm_fit %>%
forecast(h=12) %>%
autoplot() + autolayer(split_ipilsm[["test"]])
tbr_fit %>%
forecast(h=12) %>%
autoplot() + autolayer(split_tbr[["test"]])
m2_fit %>%
forecast(h=12) %>%
autoplot() + autolayer(split_m2[["test"]])
exports_fit %>%
forecast(h=12) %>%
autoplot() + autolayer(split_exports[["test"]])
imports_fit %>%
forecast(h=12) %>%
autoplot() + autolayer(split_imports[["test"]])
